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Abstract 

Sensor network attitude problem consists in retrieving the attitude of each sensor of a 
network knowing some relative orientations between pairs of sensors. The attitude of a sensor 
is its orientation in an absolute axis system. We present in this paper a method for solving 
the sensor network attitude problem using quaternion formalism which allows to apply linear 
algebra tools. The proposed algorithm solves the problem when all of the relative attitudes 
are known. A complete characterisation of the algorithm is established: spatial complexity, 
time complexity and robustness. Our algorithm is validated in simulations and with real 
experiments. 

1 Introduction 

Sensor network attitude (SNA) problem consists in retrieving the attitude of each sensors of 
a network knowing some relative attitudes between pairs of sensors. We can make the analogy 
with the sensor network location (SNL) problem which is widely studied in the literature ([2], 
[6]) and consists in retrieving sensors position from an Euclidean distance matrix eventually 
incomplete and noisy. 

Principal applications concerned are motion capture and vectorial waves measurement. In 
motion capture, the information of attitude is interesting to reconstruct the trajectory of an 
object or a body [3j. For vectorial waves, knowing the attitude allows to retrieve waves polar- 
isation. 

A basic algorithm is an algorithm which solves the SNA problem for a complete knowledge 
of the relative attitudes. As for the SNL problem [2J, having a basic algorithm allows to develop 
a more general algorithm, in particular when some relative attitudes are unknown. The focus 
of this paper is to establish a basic algorithm for the SNA problem, to characterise it and to 
validate it in simulation and in experiments. 

In section 2, we formalise the SNA problem using quaternion theory which allows to apply 
linear algebra results. 

In section 3, we propose then an algorithm solving the SNA problem for a complete and 
eventually noisy relative attitudes matrix. The most important step of the method is to esti- 
mate the highest eigenvalue and an associated eigenvector of an hermitian quaternion matrix. 
We adapt then the classical power iteration method for complex matrices [7] to hermitian 



quaternion matrices. 

In section 4, we study time and spatial complexities of the algorithm and we prove its ro- 
bustness using classical perturbation matrix results, Weyl's theorem ([5,9) and Davis-Kahan's 
theorem [1]. 

In section 5, we present an experimental validation of our algorithm using inertial systems 
placed on a polyhedron of known geometry used as a reference. 
In Appendix, we can find all mathematical missing details. 



2 Terminology for the SNA problem. 

2.1 Quaternions and Rotations. 

A quaternion q is a four components number q = go+gii+g2j+g3k where go, 9i, 92 and qs, are 
real numbers and where i, j and k are imaginary numbers satisfying i 2 = j 2 = k 2 = ijk = — 1. 
The set of quaternions is denoted by EI and can be seen as a generalisation of the complex set. 
The product defined on i, j and k induces a product on EI generalising the complex product. 
With the addition, this multiplication and the multiplication with a real, E is a non commu- 
tative real algebra and a division ring. It is important to note that a real commutes with all 
of the quaternions. 

Let q = go + gii + g2j + 53k be a quaternion, go is the real part of q. A pure quaternion is 
a quaternion with a null real part. We denote by | q| 2 := (go + 9i + 92 + oi) 1 ^ the Euclidean 
norm of q. .| 2 is a multiplicative norm. We denote by q := go — gii — g2j — g3k the conjugate 
of q. We have the following properties: pq = qp and qq = |q| 2 . 

A unitary quaternion is a quaternion of norm equal to 1. The set of unitary quaternions 
is denoted by S. A unitary quaternion q can be parametrised by an angle 9 £ [0; 2n[ and a 
unit 3D vector u := [a/3-y] T € R 3 by q = cos(0/2) + sin(#/2)u where u := «i + /3j + 7k is the 
pure and unitary quaternion deduced from u. This parametrisation allows to associate to every 
3D rotation matrix of angle and vector u, expressed in the canonical base of R 3 , a unitary 
quaternion q = cos((?/2) + sin(#/2)u using the following transformation [3]: 

9o + 9i ~ 92 - 93 2(gig 2 -g g3) 2(gig 3 + go92) 



£(q) = 



2(gig 2 + gog3) 9o-9i+92-93 2(g 2 g 3 -g gi) 
2(9193-9092) 2(g 2 g 3 + gogi) 9o-9i~92+93 



(1) 



For all unitary quaternions q, TZ(q) is then a rotation matrix and 7?.(q) = 7£(q) T = 7?;(q) _1 , 
where T is the transposition operation. This implies that the inverse of a rotation parametrised 
by q is parametrised by q _1 = q. Furthermore, for all unitary quaternions p and q, 7?.(pq) — 
7£(p)7£(q). This property shows that composing 3D rotations is equivalent to multiply the as- 
sociated quaternions. It is important to note that TZ is not an injection because 7£(q) = 7£(— q). 
Indeed, we show in appendix 7.1 that there are exactly two unitary quaternions which repre- 
sent the same rotation and they are opposed. Then, two quaternions could be far with respect 
to the Euclidean norm but they can represent the same rotation. To get over this problem, we 
can only deal with unitary quaternions with a positive real part. 

We denote by H MxJV the set of quaternion matrix of size M x N . Tr is the trace operator 
and * is the transposition-conjugation operator. An hermitian matrix is a matrix A satisfying 
A = A*. Because of the non commutativity of the algebra H, we should have to consider right 
and left eigenvalues of every matrix. We only need to consider right eigenvalues and we will not 



2 



always mention the side for further. We can show that an hermitian matrix A 6 H JVxiv has 
only real eigenvalues and is diagonalisable in an orthonormal base i.e. it exists a quaternion 
matrix U £ H JVxJV satisfying UU* = Tv such that: A = UAU* where A := diag (Ai, . . . , Ajv) 
is the diagonal matrix which contains eigenvalues of A [TO]- Finally, the Frobenius norm of 

AeH MxJV is jAj| F := Tr (^A* A) 172 - 

2.2 Attitude of sensors. 

Let 9t be the axis system of a three-components sensor and let $Ho be a reference axis system. 
The attitude of the sensor is the 3D rotation which transforms £R into 9to- Then, as a rotation 
can be represented by a unitary quaternion, the attitude of a sensor can also be represented 
by a unitary quaternion. In the sequel, we will not distinguish the attitude, the rotation and 
the unitary quaternion associated. 

We consider now N three-components sensors Si, ... , Sn with their axis system 9ti, . . . , SHjv 
and their attitude qi, . . . , qjv stored in a vector Q = [qj, . . . , qjv] T £ called the attitude 
vector. A sensor with a known attitude is called a reference. We denote by Q the sub-vector 
of Q containing quaternions associated to the references, it is the reference vector. The relative 
attitude between sensor i and sensor j is the rotation which transforms 9^ into $Hj. This 
rotation is represented by the quaternion qy = q^qj. We store those quaternions in a matrix 
denoted by O and called the relative attitudes matrix. 

2.3 SNA problem. 

Using notions and notations defined above, the SNA problem can now be enunciated as: 

"How can we estimate the attitude vector Q with an incomplete and noisy relative attitudes 
matrix O and a reference vector Q ?" 

As we precise in the introduction, we only deal with a complete relative attitudes matrix in 
this paper. In that case, the SNA problem can be seen as an inverse problem where the direct 
problem is to retrieve the relation attitude matrix from the attitude vector. This is can easily 
be done using the following formula: 

Q = QQ* (2) 
3 Basic algorithm for the SNA problem. 

We describe in this section a basic algorithm for the SNA problem based on the quaternion 
theory. First, we solve the SNA problem for a non noisy case. Then, we adapted the method 
to take into account uncertainties on numerical computations, relative attitudes and reference 
attitudes. We study the performances of our algorithm in the last subsection. 

3.1 Complete and non noisy case. 

We start by proving that every vector R £ satisfying R R* = O is unitary-right-collinear 
to the attitude vector Q. More precisely, we have the following theorem: 

Theorem 1. VR £ B N , Q Q* = RR* <=>■ 3s £ S, Q = Rs. 
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Proof. The proof of the "right to left" side is easy. Let R G § satisfying assumptions of theo- 
rem 1 and let s = Q 6 H. By multiplying by Q to the right of the equality Q Q* = RR* 

we obtain Q = Rs. Furthermore, this last equality implies that s is unitary because Q and R 
contain unitary quaternions. □ 

Theorem 1 shows that if a particular solution of equation |2| is known we can deduce the ori- 
entation vector with at least one reference. More precisely, each column of a particular solution 
is unitary-right-collinear to the attitude vector Q. Our algorithm proposes to take into account 
all of information contained in the relative attitudes matrix and the reference vector. 

We present now a method to compute a particular solution. O is an hermitian matrix 
and therefore it has real eigenvalues and an orthonormal right-eigenbase [10]. Using equation 
|2j| we can note that O is of rank 1 in EI and TrO = N. Those remarks prove that O has 
two eigenvalues : of order N — 1 and N of order 1. Furthermore, for every vector R £ §> N 
satisfying assumptions of theorem 1, we have O R = iVR R is thus an eigenvector of matrix 
O associated to the only non null eigenvalue N. Finally, estimating a particular solution of 
equation |2| is equivalent to estimating an eigenvector of O associated to the eigenvalue N i. e. 
the highest eigenvalue of O. 

We explain now how to estimate the attitude vector Q using a particular solution R and 
the reference vector Q . If there is no reference, then R can be taken as the solution but there 
is still an ambiguity due to the existence of the rotation relating R and Q. Otherwise, let R r 
be the sub-vector of R associated to the references. According to theorem 1, it exists a unitary 
quaternion s such that Q = R r s. Then s can be estimated using the following equality: 



This allows to take into account all information contained in the reference vector. The solution 
can finally be expressed as Q = (R*R r ) _1 R ( R*Q j . 



Finally, we can summarize our basic algorithm by the following list: 

step 1*) Compute an eigenvector R of O associated to N 

step 2*) Compute s 

step 3*) Compute and return Q = Rs 

This algorithm is theoretical. If we have to implement it, numerical uncertainties have to 
be taken into account. Furthermore, the relative attitudes matrix should be the issue of an 
estimation process, and would be thus an approximation of the theoretical relative attitudes 
matrix. Those remarks are also true for the reference vector. 

3.2 Complete and noisy case. 

The relative attitudes matrix and the reference vector are now considered noisy and denoted 
by Q. an d Q , respectively. We denote by Q the estimated attitude vector by the adapted al- 
gorithm described below. 

We assume that is still hermitian, with unitary quaternion elements and with diagonal 
elements equal to 1. Those assumptions are true in practice. To adapt step 1* of the theoretical 
algorithm, we have to note that N will not be, in general, an eigenvalue of O. Then, the 




(3) 
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algorithm has to find the highest eigenvalue Ai of O and an associated eigenvector R. We can 
extract R r from R and then, using Q , we can compute § as in step 2*: 

s= (r^J^^qJ (4) 

We can not directly return Q := Rs because, in general, every component will not be a unitary 
quaternion. We have to normalise each component. Finally, the adapted algorithm has the 
following form: 

step 1) Compute an eigenvector R of O associated to its highest eigenvalue Ai 
step 2) Compute s 
step 3) Compute Q 

step 4) Normalise each component of Q = Rs and return it 

4 Performances. 

We study now performances of this algorithm: spatial complexity, time complexity and 
robustness with respect to the noise. We illustrate our result in simulations under Matlab 
2008. 

4.1 Spatial and time complexities. 

First, we can note that spatial and time complexities of step 2, step 3 and step 4 are 
clearly negligible in front of those of step 1. To compute R we apply an adapted power iteration 
algorithm for hermitian quaternion matrices. The algorithm and the proof of its convergence 
are in appendix 7.3. The spatial complexity is 0(N 2 ) where N is the number of sensors in 
the network and where O is the classical Landau notation for dominated functions. The time 
complexity is experimentally estimated to be 0(N 3 ). 

4.2 Criteria and errors associated. 
4.2.1 Criteria. 

Sensor network algorithm can be seen as an optimisation problem. Indeed, it can be for- 
mulated as a minimisation of the criterion: 

Ci(P) := ||Q-pp*||f, Pes""', E r = Q r (5) 

where P r is the sub-vector of P associated to reference sensors. This optimisation problem has 
been solved for the non noisy case in section 2. As it is difficult to solve this problem in the 
noisy case, we have relaxed it in section 3 by considering two steps. First, the algorithm try 
to minimise the criterion: 

Ci(P) = ||£-PP*|||, P e H* xl (6) 

The solution retained is P := R defined as an eigenvector of O associated to its highest 
eigenvalue. Then, the algorithm minimises: 

e 2 (t):=||Q r -R r t||!, tee (7) 

The solution retained is t := s defined by equation Q. 
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4.2.2 Robustness and error bounds for step 1. 

We derive now error bounds for variables and criteria appearing in the algorithm. For every 
variable X we denote by: 

the relative error obtained for the estimation of X by X. 



In step 1 we assume that the highest eigenvalue Ai corresponds to the theoretical eigenvalue 
N. This assumption is true for low level of noise. Indeed, let Ai, . . . , Ajv be the eigenvalues 
of O indexed in the same order as eigenvalues of O as N, 0, . . . , 0. To ensure that Ai is the 
"good" eigenvalue, we have to prove that the eigenvalue associated to N can not be confused 
with others i.e. the eigen-gap d — Ai — A2 has to be strictly positive. Using Weyl's theorem, 
we show in appendix 7.2 the following inequalities: 

l-e®<^<l;Vi6[|2,W|]M<e(g) (9) 

Those properties imply that the robustness is ensured if e(0) < 1/2. In that case d > and 
Davis-Kahan theorem, adapted to hermitian quaternion matrix in appendix 7.2, leads to: 

e(O) 

This formula ensure the robustness of the step 1. 

We compute now the value of criterion Ci at point R. Let V 1 , . . . , be an orthonormal 
eigenbase for O where for all i — 1 ... N, is associated to the eigenvalue Ai. We have: 

N 

£=1>y*v: (11) 



We assume without loss of generality that step 1 returns R = V N\_ 1 . Pre-factor v N is only 
chosen to obtain concise expressions for further. Then, using definition (JsJ) and decomposition 



(111, the value of the criterion Ci for this vector can be expressed as: 



Ci(R) = 27V 2 ( 1 - ^ ) < 2iV 2 e(£i (12) 



As criterion values are absolute and not relative, the bound derived in expression ( 12 1 depends 
also on the square of the number of sensors. 



4.2.3 Robustness for steps 2, 3 and 4. 

It is difficult to explicit error bounds for e(s) and e(Q). However, robustness of step 2 and 
step 3 are ensured because they are continuous operations. The normalisation appearing in 
step 4 impacts only the angle of the estimated quaternion. Indeed, we multiply the quaternion 
by the inverse of its norm which is real. Then, this operation does not change the direction of 
the associated rotation. Finally, our algorithm is robust with respect to noise on inputs O and 

K. 



G 



We can compute the value of criterion C2 at point s: 

c 2 (g) = |[Oj[|- IIE I^;| If (13) 

\\*±r\\F 

Then, the more R r is right-collinear to Q , the more the criterion tends to and the approxi- 
mation is accurate. 

4.2.4 Simulations. 

Inputs of the algorithm are the relative attitudes matrix O and the reference vector Q . 
We suppose that errors on references are negligible in front of errors on relative attitudes. This 
assumption is true in practice. By simulations under Matlab 2008, we observe in figure [l]the 
evolution of the output error e(Q) in function of the input error e(O) limited to [0%, 10%]. 
Error are multiplied by 100 to be expressed in percent. We compare this evolution to a linear 
evolution using a linear regression. The regression line appears in figure [l] with a slope around 
of 30%. Those results show that the output error is of the same order as the input error, 
and then validate experimentally the stability of our algorithm. We trace also in figure [2] the 
evolution of criterion value versus input error. Criterion values are divided by N 2 and multiply 
by 100 to be expressed in percent. Then, our algorithm leads to acceptable criterion values. 
Furthermore, that shows that the criterion is a computable indicator of estimation quality. 

Evolution of 003) versus efO] B,Mim ° f C 1 < Q > ,ersus e <°> 




0123456789 10 

Bp) In % 



Figure 1: Estimation error. 




01 23456789 10 
6(0) in % 



Figure 2: Criterion value versus input error. 



5 Experimental validation. 

5.1 Experimental setup. 

To validate experimentally our algorithm, we use systems developed by the CEA-LETI 
called Star Watch represented in figure [3] Those systems contain two sensors, a 3-components 
accelerometer and a 3-components magnetometer, a battery and a wireless communication 
module. Analogical data coming from the sensors are sampled at 200 Hz, quantified on 12 bits 
by the system it-self and sent to an acquisition system. 
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Figure 3: A Star Watch. 
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Figure 4: Left to right. Instrumented rhombicuboctahedron by 9 Star Watch. Measures of one 
component of the accelerometer and the magnetometer of sensor 4. 

9 Star Watch are disposed on a rhombicuboctahedron represented in figure [4] We have 
registered measures of the 9 static sensors during 5 seconds at 200 Hz. One component of the 
accelerometer and one component of the magnetometer of sensor 4 are represented in figure [4] 



5.2 Relative attitude estimation. 

For each static sensor i, returned data are samples containing six measures: three com- 
ponents of the gravity field [g]pt 4 and three components of the magnetic terrestrial field [KR 
both expressed in the axis system of the sensor £Hi , where for all vector u G R 3 , [uR is the 
vector in R 3 which contains coordinates of vector u expressed in the axis system 9h. Let D\j 
be the axis system of sensor j. Then, we have: 



S(qij)[gR 
£(<to)l£k 



(14) 
(15) 



where TZ is defined in Q and qtj is the unitary quaternion associated to the relative attitude 
between sensors i and j. It is well-known in motion capture, that given equations ( 14 1 and 1 15 1, 
where [gR, [g]Kj, [hR and [h]s»tj are known, we can estimate qij [3]. For this estimation, we 
use a classical algorithm, that we call SVDQ, described in [3] and based on the transformation 
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of the equations system (14 1, (15 1 into a linear equations system. This new system is solved 
using a singular value decomposition algorithm. 



5.3 Estimation of sensors attitude. 

The absolute axis system is the axis system of the upper face of the rhombicuboctahedron. 
Sensors attitudes are known from the geometry of the rhombicuboctahedron, we store them 
in the theoretical attitudes vector Q. From Q, we compute the theoretical relative attitudes 
matrix O. It will be useful to compute the input error. 

For all measurements, which are stationary, we just conserve the average of the 5 seconds 
measurement. Then, from those measures and using SVDQ, we estimate the 36 relative atti- 
tudes between all pairs of sensors and we stored them in the estimated relative attitudes matrix 
O. Using our algorithm, we compute an estimation of the attitude vector Q where sensor 1 is 
the unique reference. 

According to section 4, we compute the input error e(O) = 2 x 10 _4 % and the output 

error e(Q) = 1.8 x 10~ 4 %. The value of the global criterion is Ci(Q) = 1.7 x 10" 5 %. Those 
results prove the applicability of our algorithm in practical situations. For information, we 
trace the errors for each pair of sensors and for each sensor in figure [5] and figure [6] 



Input errors 



(3,3) 




Output errors 



i 



T ? T 



ID 15 20 25 30 
Number of the pair (ij) (=9*(i-1)+j) 



4 5 6 

Number of the seneor 



Figure 5: Input errors. Figure 6: Output errors. 

We can observe in figure [6] that the highest output error is on sensor 6. This is consistent 
with figure [5] where highest input errors concern pair of sensors containing this sensor. This 
indicates that the algorithm conserves errors repartition. 



6 Conclusions. 

We presented an algorithm for solving sensors network attitude problem for a complete and 
eventually noisy relative attitudes matrix. The most important step relies on the estimation 
of the highest eigenvalue and an eigenvector associated of a quaternion hermitian matrix. We 
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adapted then the power iteration algorithm for this type of matrix. 

Performances of the algorithm have been studied. Spatial and time complexities are 0(N 2 ) 
and 0(JV 3 ), respectively. Robustness has been theoretically proved and validated in simula- 
tions. Simulations also give an estimation of the evolution of the output error versus the input 
error. 

Finally, the algorithm has been tested in a practical situation using attitude control system 
(Star Watch). Those results confirm the applicability of the algorithm in real situations. 

In perspective, it will be interesting to prove the efficiency of our estimator by computing 
its Cramer-Rao bound when the density of probability of input noises is known. Furthermore, 
in order to apply this algorithm to distributed sensor networks, an important perspective is to 
derive a distributed version of our algorithm where the relative attitudes matrix is incomplete. 

7 Appendix. 

In order to alleviate notations, we do not underline more variables when they are matrices 
or vectors. 

7.1 Quaternions representing the same rotations. 

To prove that it exists only two unitary quaternion which represents the same rotation, we 
derive the following theorem which links the error on quaternions and the error on rotation 
matrices associated: 

Theorem 2. V(q, q) G S 2 , \\TZ(q) - 7£(q)|| 2 = |q-q| 2v /4- |q-q|§, where |.|| 2 is Euclidean 
matricial norm. 

Proof. We have \\lZ(q) — 7£(q)|| 2 — \\I-a — 7?.(p)|| 2 where I3 is the identity matrix of order 
3 and p = qq. p is a unitary quaternion, we note r\ and v the angle and the pure unitary 
quaternion associated to p. We recall that the spectrum of a 3D rotation matrix is composed 
by l,e lr, ,e v where rj is the angle of the rotation. As I3 — 1Z(p) is a normal matrix, its 
Euclidean norm is equal to its spectral radius max {0, |1 — e"'|, |1 — e - " 7 ]} = 2\ sin(j)/2)|. We 
compute now |q — q| 2 = |1 — p| 2 . As p = cos(?7/2) +sin(?7/2)v, then |1 — p| 2 = 2| sin(r;/4)| and 
V4-|l-p|| =2|cob(jj/4)|. Finally, |q - q| 2V /4 - |q - q|jj = 2| sin(»j/2)| = ||ft(q) - TZ(q)\\ 2 . 
□ 

According to theorem 2, its easy to show that a rotation matrix is associated only to two 
quaternions which are opposed. 

7.2 Perturbation of quaternion matrices. 

In order to derive properties on eigenvalues of a quaternion matrix, we use a transformation 
which converts a quaternion matrix into a complex matrix. Then, we can apply classical results 
on perturbation matrix theory. All mathematical details on quaternion matrix theory used in 
this appendix can be found in 10 . 

A quaternion matrix M can be decomposed as M = Mi + jM 2 where Mi and M 2 are 
complex matrices in C NxN . This decomposition allows to associate to every matrix M € H JVxJV 
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a complex matrix x (M) G <£ 2Nx2N defined by: 

X is R-linear, and, satisfies |x(M)||f = \/2||M||f and ||x(M)j|2 = ||M||2 for all matrix M. 
X(M) is hermitian if and only if M is hermitian and, in that case, spectrum of x(M) is equal 
to the spectrum of M where multiplicity of each eigenvalue are doubled. Furthermore, if 
V = Vi + jVa is an eigenvector of M associated to an eigenvalue A then x(V) := [Vi V2] T is 
an eigenvector of x(M) associated to A. 

We recall now two fundamental results of perturbation matrix theory: Weyl's theorem ([9], 
0) and Davis-Kahan theorem [T] for vectorial lines: 

Theorem 3. Let A G £ N * N an d A = A + SA G £ NxN two hermitian complex matrices of 
eigenvalues «i > \a 2 \ > ... > |ojv| and |ai| > |«2 1 > •■• > |ojv|- Let d = a\ ~ &2 be the 
eigengap of 8A. We denoted by V and V unitary eigenvectors of A and A associated to oi 
and Si, respectively. Then, we have: 

Weyl's theorem. Vfc G [|1, JV|], \&k - a,k\ < \\SA\\ 2 

Davis-Kahan theorem. If d > then |sin(V,V) < ^^^ F where |sin(V,V)| : = 
^l|V-V|| F . 



As O and O are quaternion hermitian matrices, x(O) an d x(O) are hermitian matrices. 
Eigenvalues of x(O) are: N of order 2 and of order 2(7V — 1). Eigenvalues of x(O) are: 
Ai, . . . , Ajv where each eigenvalue is of order 2. x(O) and x(O) satisfy assumptions of Weyl's 
theorem which leads to: 

|Ai-AT| < ||6-0|| 2 ; ViG [|2,JV|],|Ai| < ||6-0|| 2 (17) 



By diving each inequalities in (171 by ||0||f = N and by noting that ||.||2 < ||-||f, we easily 
obtain inequalities in Inequality |Ai| < N is justified by the fact that for every matrix, its 
spectral radius is lower or equal to its Euclidean norm. 

Inequalities in Q, imply that if e(O) < 1/2 we have d > N(l — 2e(0)) > 0. Furthermore, 
we can verify that sin(x(R,), x(R-)) = e(R). With those considerations, Davis-Kahan theorem 



applied to x(O) and x(0), which satisfy the assumptions, leads to inequality (10 1 



7.3 Power iteration algorithm for a hermitian matrix. 

In the complex context, power iteration method consists in estimating the highest eigen- 
value A of a complex matrix A [?]• A has to verify |A| > \fj,\ for every other eigenvalue of 
A. The algorithm gives also an eigenvector V associated to A. We extend this theorem to 
hermitian quaternion matrix: 

Theorem 4. Let A G H JVxJV be an hermitian matrix of eigenvalues |Ai| > A2 > ... > |Ajv|- 
Let Ui, . . . , Ujv be a right-eigenbase of M N associated to A . Let Xo = J2iLi U^a; be a vector 
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in M N with ai 7^ 0. The sequence X^+i = AXi satisfies ||Xfc+i 1 1/| Xfc | — > Xi and Xfc tends 
to an eigenvector of A associated to Ai when k — > +00. 

Proof. For all integers k, X fc = A fe X = J2iLi UiA^a;, and then Ar fe X fc — > Uiai when 
k — > +00. It is easy to conclude using this limit. □ 

It is interesting to notice that the power iteration method can be extended to hermitian 
quaternion matrices because they are diagonalisable matrices and their eigenvalues are real 
and then commute with every quaternion. Furthermore, even if Ui is unknown, the condition 
ai 7^ is always satisfied due to numerical uncertainties. 

In general, to improve numerical convergence [7], the sequences considered are Vfc = 
Xfc/||Xfc|| and X^+i = A which satisfy ||Xfc+i||/||Xfc|| — > Ai and tends to an eigenvec- 
tor of A associated to Ai when k — > +00. 
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